HYRPICS: Hybrid Regularized Reconstruction for combined Parallel Imaging and Compressive Sensing in MRI

Both parallel Magnetic Resonance Imaging (pMRI) and Compressed Sensing (CS) are emerging techniques to accelerate conventional MRI by reducing the number of acquired data in the k-space. So far, first attempts to combine sensitivity encoding (SENSE) [1] imaging in pMRI with CS have been proposed in the context of Cartesian trajectories. Here, we extend these approaches to non-Cartesian trajectories by jointly formulating the CS and SENSE recovery in a hybrid Fourier/wavelet framework and optimizing a convex but nonsmooth criterion. On anatomical MRI data, we show that HYR2PICS outperforms wavelet-based regularized SENSE reconstruction. Our results are also in agreement with the Transform Point Spread Function (TPSF) or mutual coherence criterion that measures the degree of incoherence of k-space undersampling schemes. Overview: Parallel Imaging In parallel MRI, an array of L coils is used to measure the spin density ρ into the object under investigation. dl = ΣsFSlρ + nl • ρ ∈ Cn is the full FOV image, • dl ∈ Cm the signal received by each coil l, •Sl : C n → Cn denotes a sensitivity operator, •Σs : C n → Cm represents the sampling operator, Σs = [ei1, · · · , eim] ∗, •F denotes the discrete Fourier Transform, and F ∗=F−1, • nl an additive Gaussian white noise of variance σl, the between-channel covariance is assumed diagonal Λ = diag [ σ2 1Id, · · · , σ 2 LId ] . Combining Compressed Sensing and pMRI reconstruction As pointed out in [2], MR images are usually sparse in wavelet domain. Let Ψ = [Ψ1, . . . , Ψp] ∈ Cn×p, p ≥ n design a 2D orthonormal wavelet transform with ζ the wavelet coefficients such that ρ = ∑p i=1 ζ(i)Ψi = Ψζ . Compressed Sensing and Sparsity The compressive sensing theory ensures that ρ can be recovered precisely with only few observations d(k), with k ∈ {1, . . . ,m} and m ≪ n, by computing ρ̂ = Ψζ̂ where ζ̂ ∈ arg min ζ∈Cp ‖ζ‖1 s.t. d(k) = 〈ρ,φk〉 , ∀k=1:m. Hybrid regularization The reconstruction problem can be expressed as the optimization of ζ̂ ∈ argmin ζ∈Cp [JWLS(ΦΨζ) + JS(ζ) + λAJA(∇Ψζ)] . (1) • JWLS(Φρ) = ∑L l=1 σ −1 l ‖dl−(Φρ)l‖ 2 denotes the fidelity data term, andΦ = [ΣsFS1; . . . ;ΣsFSL] is the observation operator. • JS(ζ) is an l1-like prior which promotes sparsity of the recontructed solution, • JA(∇ρ) = ∑n i=1φA (√ (∂1ρ)(i) + (∂2ρ)(i) ) , φA is a Huber function and λA ≥ 0. JA is a Total Variation like prior Primal Dual Optimization We make use of the Chambolle-Pock primal-dual method [3], which has the following properties: • ”optimal” O ( 1 k ) convergence rate for the problem under consideration, • it only requires matrix-vector multiplications, • precise enough solutions in around 50 low-cost iterations. Problem (1) can be rewritten as: min x∈Cp F(Ax) + G(x). • F(y) = JWLS(y1) + λAJA(y2) with y = [y1; y2] ∈ C mL × C2n, • G(x) = JS(x) and A = [ΦΨ;∇Ψ].    yk+1 = =prox σF† { }} { (Id + σ∂F†)−1 (yk + σAx̄k+1) Dual descent xk+1 = (Id + τ∂G)−1 } {{ } =proxτG (xk − τA yk+1) Primal descent x̄k+1 = 2x k+1 − xk Correction step (due to the scheme asymmetry) Chambolle-Pock implementation. • στ = L2 where L is the highest singular value of A, • (Id + ∂F)−1(u)=argmin v∈Rn F(v) + ‖v − u‖2/2 is the resolvent (or proximal operator) of F at point u. Optimizing the undersampling scheme μ = μ(ΦΨ) = max 1≤i,j≤p,i 6=j |〈ΦΨei,ΦΨej〉| ‖ΦΨei‖ · ‖ΦΨej‖ Mutual coherence. One of the typical results relating coherence to sparsity states that if ρ = Ψζ is s-sparse and s ≤ 1+1/μ 2 then the exact recovery of ρ can be achieved by solving ζ̂ ∈ argmin ζ∈Cp ‖ζ‖1 s.t. d(k) = 〈ρ,φk〉 , ∀k=1:m. [4, 5]. The smaller the coherence, the fewer samples are needed for perfect reconstruction (noise-free case). (a) (b) (c)